package fem2;

import inf.math.LinAlg;
import math2.IsoparametricBlockGeometry;
import fem2.enu.MeshPartType;

/**
 * 
 * @author hbui
 * 
 */
public class Block extends MeshPart {

	private Node[] nodes;

	/**
	 * constructor for isoparametric block
	 * 
	 * @param nodes
	 */
	public Block(Node... nodes) {
		this.nodes = nodes;
		if (nodes.length == 4) {
			setType(MeshPartType.TET4);
		} else if (nodes.length == 8) {
			setType(MeshPartType.HEX8);
		} else if (nodes.length == 10) {
			setType(MeshPartType.TET10);
		} else {
			throw new Error("this geometry is not supported");
		}
		setBasis(BasisManager.getInstance().getBasis(getType()));
		setGeometry(IsoparametricBlockGeometry.getInstance());
		/*
		 * chec the positiveness of jacobian to ensure correct mapping
		 */
		boolean posJ = checkJacobian();
		if (!posJ) {
			throw new Error("negative jacobian is not allowed, check block");
		}
	}

	@Override
	public Node[] getNodes() {
		return nodes;
	}

	/**
	 * Computes gradients of shape functions w.r.t physical coordinates in an
	 * array
	 * 
	 * <pre>
	 * | dN1/dX1  dN1/dX2 dN1/dX3 |
	 * | dN2/dX1  dN2/dX2 dN2/dX3 |
	 * |    .        .    |
	 * |    .        .    |
	 * |    .        .    |.
	 * </pre>
	 * 
	 * @param xi1
	 *            first natural parameter
	 * @param xi2
	 *            second natural parameter
	 * @param xi3
	 *            third natural parameter
	 * @return matrix of gradients
	 */
	@Override
	public double[][] gGradientsAt(double... xi) {
		double[][] iJac = this.jacXAt(xi);
		math2.LinAlgUtilities.invert3x3(iJac);
		int n = this.getBasis().getDimension();
		double[][] gGradient = new double[n][3];
		for (int i = 0; i < n; i++) {
			double[] tmp = this.getBasis().getBasisFunction(i).gradientAt(xi[0], xi[1], xi[2]);
			double[] dNi = new double[3];
			LinAlg.multiply(false, 3, 3, 1.0, iJac, tmp, 0.0, dNi);
			gGradient[i][0] = dNi[0];
			gGradient[i][1] = dNi[1];
			gGradient[i][2] = dNi[2];
		}
		return gGradient;
	}

	public double getVolume() {
		return computeVolume();
	}
}
